home *** CD-ROM | disk | FTP | other *** search
/ Die Ultimative Software-P…i Collection 1996 & 1997 / Die Ultimative Software-Pakete CD-ROM fur Atari Collection 1996 & 1997.iso / g / gnu_c / pmlsrc23.zoo / pmlsrc / cabs.c < prev    next >
Encoding:
C/C++ Source or Header  |  1994-03-19  |  4.8 KB  |  248 lines
  1. /************************************************************************
  2.  *                                    *
  3.  *                N O T I C E                *
  4.  *                                    *
  5.  *            Copyright Abandoned, 1987, Fred Fish        *
  6.  *                                    *
  7.  *    This previously copyrighted work has been placed into the    *
  8.  *    public domain by the author (Fred Fish) and may be freely used    *
  9.  *    for any purpose, private or commercial.  I would appreciate    *
  10.  *    it, as a courtesy, if this notice is left in all copies and    *
  11.  *    derivative works.  Thank you, and enjoy...            *
  12.  *                                    *
  13.  *    The author makes no warranty of any kind with respect to this    *
  14.  *    product and explicitly disclaims any implied warranties of    *
  15.  *    merchantability or fitness for any particular purpose.        *
  16.  *                                    *
  17.  ************************************************************************
  18.  */
  19.  
  20.  
  21. /*
  22.  *  FUNCTION
  23.  *
  24.  *    cabs   double precision complex absolute value
  25.  *
  26.  *  KEY WORDS
  27.  *
  28.  *    cabs
  29.  *    complex functions
  30.  *    machine independent routines
  31.  *    math libraries
  32.  *
  33.  *  DESCRIPTION
  34.  *
  35.  *    Computes double precision absolute value of a double
  36.  *    precision complex argument, where "absolute value"
  37.  *    is taken to mean magnitude.  The result replaces the
  38.  *    argument.
  39.  *
  40.  *  USAGE
  41.  *
  42.  *    double cabs (z)
  43.  *    COMPLEX z;
  44.  *
  45.  *  PROGRAMMER
  46.  *
  47.  *    Fred Fish
  48.  *    Tempe, Az
  49.  *
  50.  *  INTERNALS
  51.  *
  52.  *    Computes cabs (z) where z = (x) + j(y) from:
  53.  *
  54.  *        cabs (z) = sqrt (x*x + y*y)
  55.  *
  56.  */
  57.  
  58. #if !defined (__M68881__) && !defined (sfp004)
  59.  
  60. #include <stdio.h>
  61. #include <math.h>
  62. #include "pml.h"
  63.  
  64.  
  65. double cabs (z)
  66. COMPLEX z;
  67. {
  68.     return( sqrt ((z.real * z.real) + (z.imag * z.imag)) );
  69. }
  70.  
  71. #else /* defined (__M68881__) || defined (sfp004) */
  72.     __asm(".text; .even");
  73.  
  74. # ifdef    ERROR_CHECK
  75.  
  76.     __asm("
  77.  
  78. _Overflow:
  79.     .ascii \"OVERFLOW\\0\"
  80. _Domain:
  81.     .ascii \"DOMAIN\\0\"
  82. _Error_String:
  83.     .ascii \"cos: %s error\\n\\0\"
  84. .even
  85. | m.ritzert 7.12.1991
  86. | ritzert@dfg.dbp.de
  87. |
  88. |    /* NAN  = {7fffffff,ffffffff}        */
  89. |    /* +Inf = {7ff00000,00000000}        */
  90. |    /* -Inf = {fff00000,00000000}        */
  91. |    /* MAX_D= {7fee42d1,30773b76}        */
  92. |    /* MIN_D= {ffee42d1,30773b76}        */
  93.  
  94. .even
  95. double_max:
  96.     .long    0x7fee42d1
  97.     .long    0x30273b76
  98. double_min:
  99.     .long    0xffee42d1
  100.     .long    0x30273b76
  101. NaN:
  102.     .long    0x7fffffff
  103.     .long    0xffffffff
  104. p_Inf:
  105.     .long    0x7ff00000
  106.     .long    0x00000000
  107. m_Inf:
  108.     .long    0xfff00000
  109.     .long    0x00000000
  110. ");
  111. # endif    ERROR_CHECK
  112. #endif    /* __M68881, sfp004    */
  113.  
  114. #ifdef    __M68881__
  115. __asm("
  116. .even
  117. .globl _cabs
  118. _cabs:
  119.     fmoved    a7@(4),fp0    |
  120.     fmulx    fp0,fp0        | x**2
  121.     fmoved    a7@(12),fp1    |
  122.     fmulx    fp1,fp1        | y**2
  123.     faddx    fp1,fp0        |
  124.     fsqrtx    fp0,fp0        | sqrt( x**2 + y**2 )
  125.     fmoved    fp0,a7@-    |
  126.     moveml    a7@+,d0-d1    | return arg
  127. ");
  128. #endif    __M68881__
  129.  
  130. #ifdef    sfp004
  131. __asm("
  132. | pml compatible lib for the atari sfp004
  133. |
  134. | Michael Ritzert, Oktober 1990
  135. | ritzert@dfg.dbp.de
  136. |
  137. | FUNCTION:    CABS(COMPLEX X)
  138. |
  139. | base =    0xfffa50
  140. |      the fpu addresses are taken relativ to 'base':
  141. |
  142. | waiting loop ...
  143. |
  144. | wait:
  145. | ww:    cmpiw    #0x8900,a1@(resp)
  146. |     beq    ww
  147. | is coded directly by
  148. |    .long    0x0c688900, 0xfff067f8
  149. | and
  150. | www:    tst.b    a1@(resp)
  151. |    bmi.b    www
  152. | is coded by
  153. |    word    0x4a68,0xfff0,0x6bfa        | test
  154.  
  155. comm =     -6
  156. resp =    -16
  157. zahl =      0
  158.  
  159. .globl _cabs
  160. .text
  161. .even
  162. _cabs:
  163.     lea    0xfffa50,a0
  164.  
  165.     movew    #0x5400,a0@(comm)    | load fp0
  166.     .long    0x0c688900, 0xfff067f8
  167.     movel    a7@(4),a0@        | load arg_hi
  168.     movel    a7@(8),a0@        | load arg_low
  169.  
  170.     movew    #0x5480,a0@(comm)    | load fp1
  171.     .long    0x0c688900, 0xfff067f8
  172.     movel    a7@(12),a0@        | load arg_hi
  173.     movel    a7@(16),a0@        | load arg_low
  174.  
  175.     movew    #0x0023,a0@(comm)
  176.     .word    0x4a68,0xfff0,0x6bfa    | test
  177.  
  178.     movew    #0x04a3,a0@(comm)
  179.     .word    0x4a68,0xfff0,0x6bfa    | test
  180.  
  181.     movew    #0x0422,a0@(comm)    | fp0 = fp0 + fp1    
  182.     .word    0x4a68,0xfff0,0x6bfa    | test
  183.  
  184.     movew    #0x0004,a0@(comm)    | sqrt(fp0)
  185.     .word    0x4a68,0xfff0,0x6bfa    | test
  186.  
  187.     movew    #0x7400,a0@(comm)    | result to d0/d1
  188.     .long    0x0c688900, 0xfff067f8
  189.     movel    a0@(zahl),d0
  190.     movel    a0@(zahl),d1
  191. ");    /* end asm    */
  192. #endif    sfp004    
  193.  
  194. #if defined (__M68881__) || defined (sfp004)
  195. # ifdef ERROR_CHECK
  196. __asm("
  197.     lea double_max,a0    |
  198.     swap    d0        | exponent into lower word
  199.     cmpw    a0@(16),d0    | == NaN ?
  200.     beq    error_nan    |
  201.     cmpw    a0@(24),d0    | == + Infinity ?
  202.     beq    error_plus    |
  203.     swap    d0        | result ok,
  204.     rts            | restore d0
  205. ");
  206. #ifndef    __MSHORT__
  207. __asm("
  208. error_plus:
  209.     swap    d0
  210.     moveml    d0-d1,a7@-
  211.     movel    #63,_errno    | NAN => errno = EDOM
  212.     pea    _Overflow    | for printf
  213.     bra    error_exit    |
  214. error_nan:
  215.     moveml    a0@(24),d0-d1    | result = +inf
  216.     moveml    d0-d1,a7@-
  217.     movel    #62,_errno    | NAN => errno = EDOM
  218.     pea    _Domain        | for printf
  219. ");
  220. #else    __MSHORT__
  221. __asm("
  222. error_plus:
  223.     swap    d0
  224.     moveml    d0-d1,a7@-
  225.     movew    #63,_errno    | NAN => errno = EDOM
  226.     pea    _Overflow    | for printf
  227.     bra    error_exit    |
  228. error_nan:
  229.     moveml    a0@(24),d0-d1    | result = +inf
  230.     moveml    d0-d1,a7@-
  231.     movew    #62,_errno    | NAN => errno = EDOM
  232.     pea    _Domain        | for printf
  233. ");
  234. #endif    __MSHORT__
  235. __asm("
  236. error_exit:
  237.     pea    _Error_String    |
  238.     pea    __iob+52    |
  239.     jbsr    _fprintf    |
  240.     addl    #12,a7        |
  241.     moveml    a7@+,d0-d1
  242.     rts
  243. ");    /* end asm    */
  244. # else    ERROR_CHECK
  245. __asm("rts");
  246. # endif    ERROR_CHECK
  247. #endif defined (__M68881__) || defined (sfp004)
  248.